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Abstract 

The finite element method has become a preeminent simulation technique in electromagnet- 
ics. For problems involving anisotropic media and metamaterials, proper algorithms should be 
developed. It has been proved that discretizing in quadratic triangular elements may lead to an 
improved accuracy. Here we present a collection of elemental matrices evaluated analytically for 
quadratic triangular elements. They could be useful for the finite element method in advanced 
electromagnetics . 
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I. INTRODUCTION 



The finite element method is a numerical technique for obtaining approximate solutions to 
boundary-value problems of mathematical physics. The method was developed and applied 
extensively for the analysis of electromagnetic problems [U El El HJ |5]. High-order vector finite 
elements have been developed which make it possible to obtain highly accurate and efficient 
solutions of vector wave equations [2J. When applied to problems involving anisotropic 
media and metamaterials, since the most of available commercial packages cannot be applied, 
proper algorithms should be developed [6|. An improved accuracy results with discretization 
in quadratic triangular elements [2J. The evaluation of elemental matrices, although is a 
simple algebra, is rather cumbersome. Here we present a collection of elemental matrices 
evaluated analytically for quadratic triangular elements. 



II. BASIC RELATIONS 



We consider an optical medium with an arbitrary cross section f2 in the xy plane. With 
a time dependence of the form exp(jut), where u is the angular frequency, from Maxwell's 
equations the following vectorial wave function is derived 

Vx([p]Vx^)-i 2 [# = (1) 

where ko is the free-space wavenumber, [p], [q] are related to the permittivity and permeabil- 
ity tensors, and <fr denotes either the electric E or the magnetic H field. Dividing the cross 
section Q into quadratic triangular elements, as shown in FigjTj we expand the transverse 
components <f> x , <fi y and the axial components <p z in each element as [2] 
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where [<f>t] e is the transverse tangential field, [4> z ] e is the longitudinal nodal field of each 
element, [U], [V], [N] are shape function vectors, and T denotes a transpose. The shape 
function [N] has six components [N] = [NiN 2 N 3 N 4: N 5 Nq\ t which are expressed in terms of 
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FIG. 1: Quadratic triangular element. 
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The area coordinates L k (k = 1, 2, 3) are given by 
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The area of the element A e is given by 
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where Xk, fjk are the Cartesian coordinates of the corner points 1 to 3 of the triangle, and 
the subscripts k, I, m progress cyclically around the three corners. The coefficients a^, bk, c k 
are given by 



Ofc = xiy m - x m yi 
b k = yi- Vm 



(6) 



The shape functions [U] and [V] are given by [2| 

/361L3 

— h^Lj, 
-Z363L1 



2A e 



[V] 



2A 



I1C2L1 
I2C3L2 

-hc\L 2 
—I2C2L3 



(7) 



with |Z fc | = + where the subscripts k, I, m progress cyclically around the three corners 
of the triangle, but I}, includes a sign which is defined differently in [21 El E] • When vector 
shape functions such as those given above are employed to represent a vector field in a 
finite element solution of vector wave equation, it is necessary to consider evaluation of the 
associated elemental matrices. These integrals can be evaluated analytically for triangular 
elements. In case of simple triangular elements, analytical relations for elemental matrices 
can be found in literature pQ . For quadratic triangular elements the integration is generally 
more involved. Analytical relations presented in [2] are rather complex. Different forms 
of elemental matrices could be useful in order to check their correctness when they are 
implemented in a developed algorithm. In the following we present a collection of elemental 
matrices for quadratic triangular elements. 

III. ELEMENTAL MATRICES EVALUATED FOR QUADRATIC TRIANGULAR 
ELEMENTS 



We used the following integration formula for the area coordinates [3] 
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Sometimes, relations becomes simpler if we take into account that 



1,^,^ = 0,1,2,3,...) 
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61 + h + h = 0, Cl + c 2 + c 3 = 
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We begin with the most simple elemental matrix, jj e [N][N] T dxdy, which can be found also 
in [21 Eli 
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d[N) d[Nf 
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d[N] d[Nf 
dy dy 



-dxdy 
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Note that Eqs.(ll) and (12) are given also in [21 USE]) but here they are simplified by using 
Eq.Q. 
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d[N] d[Nf 
dy dx 



dxdy 



3A P 



x 



36ici 
4 


-62C1 
4 


— 63C1 
4 


&2C1 





&3C1 


— biC2 
4 


362C2 
4 


— £>3 C2 
4 


&1C2 


&3C2 





-6l C3 
4 


— 62C3 
4 


363C3 
4 





&2C 3 


&1C3 


b\c 2 







bi ci + 6 2 c 2 + &3C3 


63C1 + 61C3 


&3C2 + &2C3 





&2C3 


b 3 c 2 


61 C 3 + 63 Ci 


61C1 + b 2 c 2 + &3C3 


&ic 2 + & 2 ci 







&3C1 


&2C3 + &3C2 


6 2 Ci + 61 c 2 


feiCi + 6 2 c 2 + &3C3 



// 



ox 



[Nj^^dxdy 
dy 





26i 


-62 


-63 


-61 + 26 2 


-6 2 - &3 


-h + 2b 3 




-61 


26 2 


-63 


-62 + 2&i 


-6 2 + 26 3 


—bi — 63 


1 


-61 


-b 2 


26 3 


-61 - b 2 


-63 + 26 2 


-63 + 2&i 


30 


36i 


3b 2 


-63 


8(61 + 6 2 ) 


4(6 2 + 26 3 ) 


4(6i + 26 3 ) 




-61 


36 2 


3& 3 


4(6 2 + 2&i) 


8{h + 63) 


4(63 + 26i) 




3&i 


-62 


3& 3 


4(6i + 26 2 ) 


4(63 + 26 2 ) 


8(61 + 63) 




2ci 


-c 2 


-C3 


— ci + 2c 2 


-c 2 - c 3 


-ci + 2c 3 




-ci 


2c 2 


-C3 


-c 2 + 2ci 


-c 2 + 2c 3 


-Ci - c 3 


1 


-Ci 


-c 2 


2c 3 


-Ci - c 2 


-c 3 + 2c 2 


-c 3 + 2ci 


30 


3ci 


3c 2 


-C3 


8(ci + c 2 ) 


4(c 2 + 2c 3 ) 4( Cl + 2c 3 ) 




-Ci 


3c 2 


3c 3 


4(c 2 + 2ci) 


8(c 2 + c 3 ) 


4(c 3 + 2ci) 




3ci 


-c 2 


3c 3 


4(ci + 2c 2 ) 


4(c 3 + 2c 2 ) 


8(ci + c 3 ) 



d[N} T 
dx 



-dxdy = 



12A 



x 











Zi&2(6i + 26 2 ) 


hb 2 (b 2 + b 3 ) 


Zi&2(6i + 26 3 ) 





l 2 b 2 b 3 





/ 2 6 3 (6 2 + 260 


l 2 b 3 {b 2 + 2b 3 ) 


i 2 6 3 (6i + &s) 








Z3M3 


^61(61 + 62) 


/ 3 6i(& 3 + 26 2 ) 


Z 3 6i(6 3 + 26 1 ) 





-Z1&1&2 





-Zi6i(6 2 + 26i) 


-Zi&i(& 2 + 2& 3 ) 


-Zi&i(&i + & 3 ) 








-l 2 b 2 b 3 


-l 2 b 2 (b x + 6 2 ) 


-/ 2 6 2 (6 3 + 26 2 ) 


-Z 2 6 2 (6 3 + 26 1 ) 


-/3&1&3 








-Z 3 & 3 (&i + 26 2 ) 


-kh{b 2 + 63) 


-^3& 3 (&1 + 26 3 ) 



1 



dx 



hb\c 2 





-Z3&1C3 






-Z162C1 





12A e 



hhci 


~hhc2 



x 



ZiC2(6i + 26 2 ) Zic 2 (& 2 + 63) ilC2(6l + 26 3 ) 

l 2 c 3 (b 2 + 26i) Z 2 c 3 (6 2 + 26 3 ) Z 2 c 3 (&i + 63) 

Z 3 ci(fei + & 2 ) / 3 ci(6 3 + 26 2 ) Z 3 ci(fe 3 + 26i) 

-Zici(62 + 26i) -Z lCl (6 2 + 26 3 ) -iMh + h) 

-l 2 c 2 {bi + b 2 ) -l 2 c 2 (b 3 + 2b 2 ) -Z 2 c 2 (6 3 + 2b ± ) 

-Z 3 c 3 (6 1 + 26 2 ) -l 3 c 3 (b 2 + b 3 ) -Z 3 c 3 (&i + 26 3 ) 



[C/]3^cbdy - 



hb 2 ci 





-Z363C1 





hhc 2 


-hb\c 2 





12A e 



khc 3 




x 



Z 1 6 2 (c 1 + 2c 2 ) 
Z 2 6 3 (c 2 + 2ci) 

^3&l(ci + C 2 ) 

-Zi6i(c2 + 2ci) 
-l 2 b 2 c 3 -l 2 b 2 (ci + c 2 ) -l 2 b 2 (c 3 + 2c 2 ) 
-Z 3 6 3 ( Cl + 2c 2 ) -Z 3 6 3 (c 2 + c 3 ) 



Zi& 2 (c 2 + c 3 ) 
l 2 b 3 (c 2 + 2c 3 ) 
Z 3 6i(c 3 + 2c 2 ) 
-Z 1 6 1 (c 2 + 2c 3 ) 



Z 1 6 2 (c 1 + 2c 3 ) 

Z 2 &3(C1 + C 3 ) 

Z 3 6i(c 3 + 2 Cl ) 
-Zi6i(ci + c 3 ) 
-l 2 b 2 {c 3 + 2ci) 
-Z 3 6 3 ( Cl + 2c 3 ) 



[y]3^dad/y 



1 



ZiC!C 2 






-^cic 3 





hc 2 c 3 


-kc\c 2 





12A 



Z3C1C3 




x 



Zic 2 (ci + 2c 2 ) 
/2C 3 (c 2 + 2ci) 

Z3Cl(ci + C 2 ) 

-ZiCi(c 2 + 2ci) 
-Z 2 c 2 c 3 -Z 2 c 2 (ci + c 2 ) 



Zic 2 (c 2 + c 3 ) 
/2C 3 (c 2 + 2c 3 ) 

Z 3 Ci(c 3 + 2c 2 ) 

-ZiCi(c 2 + 2c 3 ) 
-Z 2 c 2 (c 3 + 2c 2 ) 



Zic 2 (ci + 2c 3 ) 

Z2C 3 (ci + C 3 ) 
^Cl(c 3 + 2ci) 

-ZiCi(ci + c 3 ) 
-kc 2 {c 3 + 2ci) 
-Z 3 c 3 (ci + 2c 2 ) -Z 3 c 3 (c 2 + c 3 ) -Z 3 c 3 (ci + 2c 3 ) 



JJ[U][U} T dxdy 



48A f 



x 





hhb 2 b 3 


hhhfo 


-lihh 


-hhb\ 


-2hl 3 b 2 b 3 


hhb 2 b 3 


2l\b\ 


hhhh 


-2hl 2 b x b 3 


-l\b 2 b 3 


-hhb 2 3 


hhbib 2 


kkhh 


2l\b\ 


-khb\ 


-2l 2 l 3 b x b 2 


-l 2 3 hb 3 




-2l x l 2 b x b 3 


-khb\ 


2l\b\ 


hl 2 bib 2 


hhhh 


—l\l 2 b\ 


-l\b 2 b 3 


—2l 2 l 3 b\b 2 


l\l 2 b\b 2 


2l\b\ 




-2Z1Z36263 


-hhb\ 


-l 2 3 hb 3 




hhb 2 b 3 


2l\b\ 



6 



JJ[V][V} T dxdy 



48A 



x 



2h\c\ 


^l/2C 2 C 3 


hhc\c 2 


-l\c\c 2 


-hhc\ 


-2hl 3 c 2 c 3 


( '1 ( '2 < ^2^3 


21 2 r 2 


In /o O-i On 




t 2 L 2 L 3 


—hhc 2 


hhc\C 2 


hhc\C 3 


2l\c\ 


-hhcl 


-2/ 2 / 3 CiC 2 


-l 2 3 c\c 3 


— l\c\C 2 


—2l\l 2 c\C 3 


L i b 3 ^\ 


2l\c\ 


<-l<-2 , ^l , ^2 


/1 /ofi fo 

< , l<'3 , ^l ( ^3 




-l 2 c 2 c 3 - 


-2/oZqCi Co 


/1 loCi Co 


2l 2 c 2 


v Z L 6^Z ^0 


-2hl 3 c 2 c 3 


-hhc\ 


— ZoTi Tq 
3 ^ "J 




"Z^o^Z *-<6 


2/3C3 




X 








2ljb 2 c 2 


hhhc 3 


hl 3 b 2 ci 


-l\b 2 c x 


—hhb 2 c 2 


-2hl 3 b 2 c 3 


hkhc2 


2l\b 3 c 3 


hhhci 


-2hl 2 b 3 ci 


-l\b 3 c 2 


~hhhc3 


hhhc 2 


hhhc 3 


2l\b x c x 


-hhhci 


-2l 2 hb\c 2 


-l 2 3 b\c 3 


-l\b\c 2 


-2lil 2 b\c 3 


-hhb\C\ 


2ljb lCl 


hhb\c 2 


hhb\c 3 


-hl 2 b 2 c 2 


-l\b 2 c 3 


-2l 2 l 3 b 2 ci 


hhb 2 ci 


2l 2 b 2 c 2 


hhb2C 3 


-2/1/363C2 


-hhbzc 3 


-l\b 3 ci 


hhb 3 ci 


hhbzc 2 


2l 2 3 b 3 c 3 



d\U] d[Uf 
dy dy 



-dxdy 



1QAI 



a n 


dl2 


«13 


a i4 


Ol5 


a 16 


au 


0-22 


a 2 3 


a 2 4 


025 


026 


an 


a 2 3 


a 33 


a 34 


^35 




a u 


a 24 


a 34 


a 44 


a 45 


a 46 




025 


035 


a 4 5 






a i6 


«26 


«36 


a 4 6 


«56 


^66 



an = 


;2t,2 2 

hb 2 c x , a 12 


= hhbzb 3 c x c 2 , a 43 


= hhbib 2 C!C 3 , 


a i4 = 


-l\b\b 2 c x c 2 , 


Ol5 = -/l/2&2 c l c 3> 


Oi6 = -/1/3&2&3C1, 


a 22 = 


i2i2 2 

L 2°3 C 2i °23 


= hhb\b 3 c 2 c 3 , a 24 


= -hl 2 bib 3 cl, 


^25 = 


-l\b 2 b 3 c 2 c 3 , 


«26 = -hhb\c\c 2 , 


;2r,2 2 
033 = / 3 »1 C 3> 


a 34 = 


-hhblc 2 c 3 , 


«35 = -hhb\b 2 c\, 


^36 = -/3O1O3C1C3, 


a 44 = 


i2i2 2 

L l°l C 2i °45 


= hhbib 2 c 2 c 3 , a 46 


= hhhhc^, 


^55 = 


;2j2 2 

L 2°2 C 3i °56 


= hhb2b 3 cic 3 , a m 


- l 2 b 2 r 2 

— l 3 u 3 c l- 



The matrix JJ e 8 ^ dxdy is expressed like in Eqs.(|23|) and (|24|) by applying the trans- 
formation 5<->c. One obtains 
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